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1. Introduction 



The equation of state(EoS) of QCD, i.e. thermodynamic quantities such as the pressure {p), 
energy density (e), or entropy density {s) as a function of the temperature, is important in under- 
standing the high-temperature behavior of QCD. The EoS is not only of theoretical interest, but is 
directly applicable to the dynamics of the quark-gluon plasma (QGP), whether in the context of 
interpreting the results of heavy-ion experiments or modelling the behavior of hot, dense matter in 
the early universe. 

Currently, some aspects of the EoS are accessible at high temperatures via perturbation theory 
[|l|], or at low temperatures via the hadron resonance gas model (HRG) However, a truly non- 
perturbative, first-principles calculation can only be done using lattice QCD. 

Recently, there have been several detailed studies of the QCD EoS using improved actions and 
physical (or almost-physical) values for the light and strange quark masses [|3|, ^, ^. While these 
studies accurately represent the three lightest quark flavors, they, like all previous calculations, 
neglect the effect of the charm quarks. However, recent work ^ based on perturbation theory 
indicates that the charm quark contribution may start to become significant at temperatures T ~ 350 
MeV, well within the temperature range covered by these latest EoS calculations. 

In this work, we present a partially quenched study of the charm quark contribution to the 
QCD equation of state using p4fat3 fermions on a 2-i-lf gauge background in the temperature range 
0.90rc <T < 4.2Tc. Currently, the calculations are done only for Nt =4 lattices, with the intention 
of continuing to A'^r = 6 and possibly Nf = S. 



2. Equation of State on the lattice 

We will first review the formalism for calculating the equation of state on the lattice via the 
"integral method", as described in [0]. First, consider the grand canonical partition function of 
QCD at vanishing quark chemical potential (Hq = 0): 

Z(V,T) = I [^A^][W][W] exp {- j^d\ ^qcd^ (2.1) 

We can extract various thermodynamics quantities from \nZ{T,V) such as the grand canonical 
potential {Q.{T,V)), pressure (p), or energy density (£): 

niT,v) = rinz(r,y); e = ^ = P = (2.2) 

^ ' ^ ^ ' ^ V Vd{l/T) ' dV 

In the thermodynamic limit, the grand canonical potential is an extensive quantity {Q. ~ V). Thus, 
we can write: 

p = -\xvZ{Ty)- e = — (2.3) 

We can also define another quantity, called the interaction measure /, which is the trace of the 
energy-momentum tensor: 

L = ^_Z}1 = T±{^) (2 4) 

j4 j4 dT 
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Unfortunately, calculating the pressure by computing In Z(r,y) directly is not possible using 
stochastic Monte Carlo methods. However, we are able to compute certain derivatives of lnZ(r, V) 
which, in the lattice formulation, become the expectation values of operators. In particular, we can 
obtain the gauge action and the chiral condensate in this manner: 

, , d\nZ I , dXrvZ 

&> = ^; (w.> = ^ (") 

where j8 = 6/go is related to the bare coupling go and ihq denotes a bare quark mass. 

On the lattice, finite temperatures can be simulated by limiting the temporal size of the lattice. 
This is directly analagous to restricting the integral in ^]l| to a finite interval in the imaginary time 
direction. The temperature is related to the temporal extent by T^' = Nta{gQ,mq{gQ)), where the 
lattice spacing a is a function of both go and nig, and the bare quark masses must be adjusted as 
a function of go along some renormalization group trajectory to keep physical quantities fixed, i.e. 
a line of constant physics. 



We can rewrite the expression g4| for the interaction measure in terms of lattice quantities: 



I fN^'f dp ((s,{T))-{S,m)+l^^({W,{T))-{wM)] (2-6) 



Here, the sum over the index q indicates a sum over all quark flavors. We have also normalized the 
interaction measure by subtracting the T = contribution, disentangling the vacuum contribution 
from the thermal effects. 



Using /, we can reconstruct other thermodynamic quantities. For example, we use 2.4 to 
express the pressure in terms of lattice quantities: 



P 

y4 




(2.7) 



^=(1) £(^((5,(n)-(5,(0)>)+l||f((w.(n>-('l^^^^^ (2.8) 

Once we have /(r) and /?(r), we can easily reconstruction the energy density, £(r), or the entropy 
density, ^(r). For a more detailed discussion of our EoS calculation with 2+1 flavors, see 

3. Calculation Method 

For this calculation, we are concerned only with the charm quark contribution. The relevant 



part of 2.6 is: 

To calculate ^\^c, we make partially quenched measurements of the chiral condensate on pre- 
viously generated 2+ If gauge configurations, both at finite temperature (A'; = 4) and zero tempera- 
ture (l^t = 32). These dynamical configurations use the p4fat3 fermion action , and a tree-level 
improved Symanzik gauge action. 
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In order to determine the temperature of each ensemble, we have chosen to set the scale using 
the static quark potential parameter r^. ro is defined as: 

f^dV^) ^1.65 (3.2) 

We can then convert to physical units by using tq = 0.469(7) fm. [^. Note that the bare quark 
masses (m„j and m,) are tuned as a function of j8 so that physical quantities such as nij^rQ, mrfrQ, 
and rriKro take on approximately constant values over the entire temperature range. 

To determine Ic, we not only need on zero and finite temperature lattices, but also 
dp/d\na. We can deduce dp /dlna from our scale-setting calculations of tq: 

dp _/P _( dln{a/ro) \-' 
dlna da \ dp J 

To get a smooth function for dp /dlna, we fit ^(j3) to a renormalization-group inspired ansatz. 
Table 1 gives details of the input parameters for the lattices, the number of trajectories for the 
different ensembles, as well as ro and T /T^. For more details about scale setting and the line of 
constant physics, see [Q]. 

















Traject 


Dries 




Volume 


mud 


ro/a 


T/T, 


mc{r]c) 


m,(7/»F) 


Af, = 32 


Nt=4 


3.277 


16^ 


.00765 


1.797(19) 


0.90 


231 


2.25 


3250 


12160 


3.335 


16^ 


.00570 


2.033(17) 


1.06 


1.57 


1.42 


2630 


14280 


3.351 


16^ 


.00592 


2.069(12) 


1.10 


1.50 


1.35 


6950 


12420 


3.382 


16^ 


.00520 


2.225(13) 


1.20 


1.22 


1.11 


2270 


8110 


3.41 


163 


.00412 


2.503(18) 


1.31 


1.02 


.858 


2790 


16000 


3.46 


163 


.00313 


2.890(16) 


1.50 


.670 


.650 


2510 


10200 


3.49 


163 


.00290 


3.223(31) 


1.62 


.566 


.529 


4290 


9420 


3.51 


163 


.00259 


3.423(61) 


1.70 


.508 


.473 


2450 


10000 


3.54 


163 


.00240 


3.687(34) 


1.83 


.446 


.417 


4060 


6250 


3.57 


243 


.00212 


4.009(26) 


1.98 


.386 


.347 


2460 


21190 


3.63 


243 


.00170 


4.651(41) 


2.28 


.304 


.288 


3290 


10000 


3.69 


243 


.00150 


5.201(48) 


2.61 


.257 


.244 


2290 


9470 


3.76 


243 


.00130 


6.050(61) 


3.05 


.213 


.205 


1110 


33370 


3.82 


243 


.00110 


6.752(96) 


3.46 


.181 


.190 


3000 


35000 


3.92 


243 


.00092 


7.59(12) 


4.23 


.162 


.157 


4080 


35870 



Table 1: Input parameters for the ensembles on which we have performed measurements. Note, nis = 10m„f/ 
for all of these ensembles. Also given are ro, T /Tc, nic determined from m-q^ and mjjm, and the number of 
trajectories for Ni =4 and A', ~ 32. 

4. Setting the physical charm mass 

In addition to dP /dlna, we also need drhc/dp to calculate 7^. Thus, we need some method to 
determine the bare charm quark mass so that it also sits on a line of constant physics as we vary p. 
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Figure 1: On the left, we see ihc/ihud as a function of j3. On the right, we see dihc/dji. 



Ideally, if we were simulating in the scaling regime for charm quark observables, the ratio 
fhc/thud should be a constant throughout the entire temperature range and independent of the quan- 
tity chosen to fix the. Unfortunately, rhc is quite large on the available lattices, so the ratio rhc/rhud 
depends on j8 and the particular quantity used to set . 

In order to understand the extent of this scaling violation, and how it affects our calculation, 
we have chosen to fix rhc by setting two different charmonium states to their physical masses {m-q^ 
= 2.980 GeV and trij/^, = 3.097 GeV). This produces a range of values 0.16 < mc < 2.4, with 
exact figures given in Table 1. As seen in Figure |l|, we find fairly reasonable scaling {rhc/rhud ~ 
170) for the finer lattice spacings (j3 > 3.5 — > T > l.lTc). However, as we lower j8 onto coarser 
lattices, we find that this approximate scaling breaks down, and this ratio rhc/rhud begins to increase 
dramatically. We also find that using the r\c state tends to give a systematically higher value for rhc, 
although this ultimately does not have a large affect on the final calculation. 

Once we have ihc{p), we can fit this data to a RG-inspired ansatz for the running of the bare 
quark mass. This allows us to calculate a smooth curve for dnic/d^, also shown in Figure [1} We 
see that this j8 -function increase s dramatically as we move to smaller j8, in contrast to the scaled 
version of drhud/d^, which shows a much milder j8 -dependence. 

5. Results 

Figure § shows our results for the charm contribution using both r\c and J /^' to determine rhc, 
as well as a free-field calculation using the physical value of the charm mass nic ~ l.lGeV. The 
choice of rjc or does have some effect, but both curves share the same qualitative features. 
Notably, the interaction measure Ic increases drastically as we move through the transition. This is 
consistent with the notion that the charmonium states are the lightest charm states that contribute 
below the transition, and thus are exponentially suppressed by their heavy mass. 

It is interesting to note that the value of Ic already becomes quite large even just above the 
transition (T ~ 1.1 T^), decaying away slowly as the temperature is increased. This is in contrast to 
the free-field calculation, where the interaction measure peaks at T ~ 0.35mc ~ 400MeV, although 
of course the free-field result knows nothing about the crossover transition. Perhaps more surprising 
is the fact that Ic is significantly greater than the free-field value for T < 3Tc. As a result, p and e 
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(e-3p)/T^ (e-3p)/T' 




Figure 2: In the left column, the charm contribution to the interaction measure 

pressure (^), and 

energy density (^), respectively. Diamonds use r\c while the triangles use J to set ntc- The solid curve 
gives the free-field results, while the solid bars on the upper right give the T ^ °o limit. In the right column, 
the charm contribution (stars) is added onto the pure 2+lf part (diamonds). 

increase much faster just above the transition than expected. This causes p and e to overshoot their 
continuum Stefan-Boltzmann values T — > oo by a large amount, even at finite temperature. This 
difference may be attributed to finite lattice spacing corrections at A^f = 4. Indeed, for the 2+lf 
calculation, the contribution of the quark condensates to {e — 3p)/T^ decreases going from Nt =4 
to Nt = 6. 

Figure ^ also shows the charm contribution superimposed on the 2+ If calculation at A^, = 4. 
Although we see little change in the interaction measure because of the high 2+ If peak, the pressure 
and energy change by quite a large amount for T > T^. 

6. Conclusion 

We have made a partially quenched calculation of the charm quark equation of state. Current 
EoS studies are in a temperature regime that is a significant fraction of the charm quark mass. 
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where the dynamics of the charm quark may have a significant effect and can no longer be ignored. 

Admittedly, this calculation has unquantified, possibly large systematic errors. First of all, 
our calculations are partially quenched, so that the full dynamics of the charm quark may not be 
accurately reflected in their "back-reaction" on the gauge fields. To the extent that the charm quark 
is heavy, this should be a small effect, but may be problematic at higher temperatures and at finer 
lattice spacings. Secondly, rhc is a large fraction of, or even exceeds, the lattice spacing on many 
of the lattices on which we measure. This introduces large cut-off effects, and it is appropriate to 
ask whether the p4 fermion formulation can give even a reasonable facsimile of the charm quark 
dynamics at these heavy masses. 

As there are no plans for large-scale simulations involving a dynamical charm quark, it seems 
we must tolerate the first problem. The second problem may be alleviated somewhat be moving to 
finer lattice spacings, {Nt = 6, 8). This may give us a better grasp of the discretization errors in the 



current calculation. Furthermore, the recent development of the HISQ action ||10p, which, among 
other improvements, removes the leading order 0{ma)^ corrections to the dispersion relation at 
tree-level, provides another possible tool to mitigate the large cut-off effects currently present. 

The future plans for this calculation is extension to Nt = 6 and possibly A/, = 8, as well as 
investigating the possibility of adapting the methods used by the HISQ action. Hopefully, this will 
allow us to better understand why the Nt = 4 calculation so badly overshoots the Stefan-Boltzmann 
limit, and to better quantify the cut-off effects. 
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